function [cust_R] = custom_reaction(xx,t)
    %return cust_R with the value of the reaction coefficient "R" for each node
    [nnod,ndm] = size(xx);
    cust_R = zeros(nnod,1);
    %% completar desde aqui
end